# -*- coding: utf-8 -*-
"""This module contains numeric manipulation functions which don't fit within
any of the other modules in the numbers package.
"""

"""Project Euler Solutions Library

Copyright (c) 2011 by Robert Vella - robert.r.h.vella@gmail.com

Permission is hereby granted, free of charge, to any person obtaining a copy
of this software and associated documentation files (the "Software"), to deal
in the Software without restriction, including without limitation the rights
to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
copies of the Software, and to permit persons to whom the Software is
furnished to do so, subject to the following conditions:

The above copyright notice and this permission notice shall be included in
all copies or substantial portions of the Software.

THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
THE SOFTWARE.
"""

import math
from itertools import takewhile
from euler.core import infinity, LazyEvaluatedList
from euler.numbers.primes import Primes, prime_factors
from euler.list_functions import join, true_product, key_max_value_pairs, \
                                key_count_pairs

def factors(number, primes = None):
    """Returns all the factors of [number], including 1 and [number] itself.
        
    Optional Parameters:
        primes - An already active prime cache 
                    (See: euler.numbers.primes.Primes), as this algorithm 
               depends on the generation of prime numbers, passing an instance
               of Primes, with a number of pre-generated results, will actually
               reduce execution time. New primes generated within this function
               will also be added to the cache.
    """              
    
    #Initialise cache if it has not been specified.
    primes = primes or Primes()
    
    #The prime factors of [number].
    prime_factors_of_number = set(prime_factors(number, primes))
    
    #The list of factors.
    factors = []
    
    #Generate, and add, all the multiples of the unique prime factors of 
    #[number], which are less than or equal to its square root.
    number_square_root = int(math.sqrt(number))
    
    for prime_factor in prime_factors_of_number:
        factors += filter(lambda x : number % x == 0,
                              takewhile(lambda m : m <= number_square_root,
                                            multiples(prime_factor)
                                        )
                              )
                       
    #Add factors above square root, by taking the quotients of [number] and
    #all the factors which have been generated so far.
    factors += map(lambda x : int(number / x), factors)
   
    #Add 1 and [number] to the factors.
    factors += [1, number]
   
    #Return unique factors
    return list(set(factors))
    
        
def number_of_divisors(n, primes = Primes()):
    """Returns the number of factors of [n]. This function is a lot faster than
    counting the result of the function, factors.    
    
    Optional Parameters:
        primes - An already active prime cache 
                    (See: euler.numbers.primes.Primes), as this algorithm 
               depends on the generation of prime numbers, passing an instance
               of Primes, with a number of pre-generated results, will actually
               reduce execution time. New primes generated within this function
               will also be added to the cache.
    """
    return true_product(instances_of_prime + 1 
                         for prime, instances_of_prime in 
                            key_count_pairs(prime_factors(n, primes)))

  
def factorial(number):
    """Returns [number]th factorial, starting with 1 at 1."""
    return true_product(range(number, 0, -1))

def multiples(number):
    """Generates all the multiples of [number]."""
    return infinity(number, number)

def lcm(numbers, primes = Primes()):
    """Returns the least common multiple of the numbers in the iterable, 
    [numbers].
    
    Optional Parameters:
        primes - An already active prime cache 
                    (See: euler.numbers.primes.Primes), as this algorithm 
               depends on the generation of prime numbers, passing an instance
               of Primes, with a number of pre-generated results, will actually
               reduce execution time. New primes generated within this function
               will also be added to the cache.
    """
    
    #Calculate the grouped prime factors of [numbers]. That is, generate the
    #prime factors of each number in the form (prime, instances of prime).
    #
    #Example:
    #
    #numbers = [4, 30, 10] : grouped_factors = ((2, 2), (5, 1), (2, 1), (3, 1), 
    #                                            (5, 1), (2, 1))
    all_factors = (prime_factors(n, primes) for n in numbers)

    
    grouped_factors = (key_count_pairs(factors_of_number) for \
                          factors_of_number in all_factors)

    grouped_factors = join((list(factors) for factors in grouped_factors))

    #Return the product of the primes to the power of their instances. Be sure
    #to only take the the maximum number of instances for each prime.
    #
    #Example:
    #
    #grouped_factors = ((2, 2), (5, 1), (2, 1), (3, 1), 
    #                  (5, 1), (2, 1))
    #
    #returned value = 2^2*3^1*5^1 = 60    
    result = true_product(item[1][0] ** item[1][1] for item in \
                     key_max_value_pairs(grouped_factors, lambda n : n[0]))

    return result

def quadratic_roots(a, b, c):    
    """Returns the quadratic roots for the formula [a]x**2 + [b]x + [c]."""
    
    #This is a python based implementation of the standard equation for the 
    #roots of a quandratic formula.
    b_squared_minus_4ac = b ** 2 - 4 * a * c
    
    if b_squared_minus_4ac <= 0:
        return None
    
    root_b_squared_minus_4ac = b_squared_minus_4ac ** 0.5
    two_a = 2 * a
    
    return ((-b - root_b_squared_minus_4ac) / two_a, \
            (-b + root_b_squared_minus_4ac) / two_a)


def totients(maxn):
    """Generates φ(n) for all n <= [maxn]. φ(n) gives the number of positive 
    integers which are < n and coprime to n.
    """
    
    #This algorithm is based on the following assumptions, it is important to
    #take note of them before reading further, so as to better understand the
    #algorithm.
    #
    #φ(n) * φ(k) = φ(n * k)
    #φ(n) = n - 1 if n is prime
    #φ(n**k) = φ(n)*n**(k - 1)     
    
    #Yield φ(1), which is 1.
    yield 1
    
    #The cache which holds all φ(n) for n = 2 to maxn.
    cache = [0 for i in xrange(2, maxn + 1)]
    
    #The cache which holds all the primes and their exponents, as well as the 
    #product of their exponents and the exponents of previously calculated 
    #primes, including the current prime itself.
    #
    #if [maxn] = 20, this cache will contain the following values by the time 
    #the prime, 11, is calculated:
    #
    #[(2, [4, 8, 16]), (3, [6, 9, 12, 18]), 
    # (5, [10, 15]),  (7, [14]), (11, [])].
    primes_prime_exponents = []
    
    def process_exponent(exponents_list, exponent, totient):
        """Processes the exponents of a prime and adds their products with 
        the previously calculated prime exponents to [exponents_list].
        
        Parameters:
            exponents_list - An array to which the values calculated within 
                            this function will be added.
            
            exponent - The exponent for which the values will be calculated.
            totient - φ(exponent).
            
        """
        
        #Store the result of φ(exponent) within the cache.
        cache[exponent - 2] = totient
        
        #Go through all the previously calculated values within 
        #primes_prime_exponents, and add new values based on the product of the
        #current exponent with these values.
        for prime_prime_exponents in primes_prime_exponents:
            current_exponent_times_prime = exponent * prime_prime_exponents[0] 
            if current_exponent_times_prime > maxn:
                break
            
            cache[current_exponent_times_prime - 2] = totient \
                                         * cache[prime_prime_exponents[0] - 2]
            exponents_list.append(current_exponent_times_prime)
            
            
            if prime_prime_exponents[1] is not None:
                for prime_exponent in prime_prime_exponents[1]:
                    current_exponent_times_exponent = exponent * prime_exponent
                    
                    if current_exponent_times_exponent > maxn:
                        break
                                
                    cache[current_exponent_times_exponent - 2] = totient \
                                                * cache[prime_exponent - 2]
                    exponents_list.append(current_exponent_times_exponent)

    #For all 2 <= n <= maxn.
    for n in xrange(2, maxn + 1):
        #If φ(n) has not been calculated.
        if cache[n - 2] == 0:    
            #As the process calculates φ(n) for all exponents of a prime
            #and the product of those exponents with previously calculated
            #prime exponents, it follows that if φ(n) has not yet been 
            #calculated then n is prime.       
            totient = n - 1
            new_prime_exponents = []
            
            #Process the first exponent of n which is n itself.
            process_exponent(new_prime_exponents, n, totient)
            
            #Process the remaining exponents of n.
            current_exponent = n * n
            totient *= n
            
         
            while current_exponent <= maxn:
                new_prime_exponents.append(current_exponent)
                            
                process_exponent(new_prime_exponents, current_exponent, totient)
                
                current_exponent *= n
                totient *= n
            
            #If we have an exponent of n or a product of an exponent of n 
            #and a previously calculated exponent, then sort this list and
            #store it in the cache, prime_prime_exponents. Otherwise, store
            #n by itself.               
            if len(new_prime_exponents) > 0:
                new_prime_exponents.sort()
                primes_prime_exponents.append((n, new_prime_exponents))
            else:
                primes_prime_exponents.append((n, None))
        
        #Yield φ(n).   
        yield cache[n - 2]
         
        
  
    
class Factorials(LazyEvaluatedList):
    """A cached list of factorials."""
    
    def __init__(self):
        LazyEvaluatedList.__init__(self, [1, 1, 2, 6, 24, 120])
  
    #Override.
    def _next_item(self, last_result, number_of_known_terms):
        return self[-1] * len(self.terms())

